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and distribution (Pecl et al., 2017), which endangers the existence of multiple organisms 
(Chen et al., 2011; Trew and Maclean, 2021). 


Climate-driven changes in species’ potential distribution can be assessed using species 
distribution model (SDM) (Zurell et al., 2023). SDM can be used to predict suitability 
of habitats for different species under different environmental conditions making it an 
extremely significant tool for conservation planning and decision-making (Rahman et 
al., 2019). This modeling approach is especially beneficial for species that are vulnerable 
to climate change, such as the Mindanao horned frog. The predictive capability of SDM 
enables the targeted allocation of conservation areas, highlighting the significance of 
focusing conservation efforts in areas of high conservation value, especially for endemic 
species with restricted ranges (Guisan ef al., 2013). 


Amphibians, due to their distinct physiological requirements and life cycles, are 
exceptionally susceptible to the consequences of climate change, making them one of 
the most endangered animal groups (Hof et al., 2011; Li, Cohen and Rohr, 2013; Luedtke 
et al., 2023). Fluctuations in temperature, changes in rainfall patterns, and an increased 
occurrence of severe weather events greatly affect the ability of amphibians to survive, 
reproduce, and spread (Li, Cohen and Rohr, 2013; Campbell, Miller and Muths, 2020; 
Rollins-Smith and Le Sage, 2023). These impacts are intensified for endemic species 
whose distribution is restricted and limited (Manes ef al., 2021) urging the call for 
additional studies for this vulnerable group (Luedtke ef al., 2023). 


The Mindanao horned frog (Pelobatrachus stejnegeri) is an endemic amphibian species 
found exclusively on the island of Mindanao in the Philippines. It is known for its unique 
ecological role and its importance as a bio-indicator (IUCN SSC Amphibian Specialist 
Group, 2020). This terrestrial frog is commonly found in both main and secondary 
montane and lowland rainforests having mountain streams as the primary breeding 
grounds. Being endemic, these species are extremely susceptible to habitat alterations 
resulting from climate change, mostly due to their specific biological requirements and 
restricted geographic range (Bickford et al., 2010; Manes et al., 2021). However, no 
studies have been conducted yet, exploring the habitat suitability and distribution of the 
Mindanao horn frog and assessing its response under changing climates. 


This study aims to model the distribution of the endemic Mindanao horned frog (P. 
stejnegeri) under current and future climatic conditions. Additionally, it aims to measure 
the impact of climate change on this species' geographic distribution through 
determination of regions where the habitat may become more or less suitable based on 
future climate scenarios. 


Methodology 

Study Area and Occurrence Records 

The study was carried out in the island of Mindanao in the Philippines, where P. 
stejnegeri is endemic (Figure 1). Mindanao is located in the southern part of the 
Philippines and is one of the three island groups of the country, covering an approximate 
land area of 95,468.17 square kilometers (PhilAtlas, 2022). The occurrence points of P. 


stejnegeri were collected from the Global Biodiversity Information Facility 
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(https://www.gbif.org/) and survey data conducted in 2012—2023 by the author. 
Duplicate geographic coordinates and coordinates falling outside the study area 
boundary were removed from the data points to clean them up. Using SDMtoolbox 
(Brown, Bennett and French, 2017), the distribution data were subjected to spatial 
rarefaction at a 5 km distance in order to eliminate sampling bias and spatial 
autocorrelation. 
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Figure 1: Location of the study area with observations of P. stejnegeri (red dots). 
Environmental Variables 


This study utilized 21 environmental layers with a spatial resolution of 30 seconds to 
create a distribution model of P. stejnegeri. These include topographical factor 
(elevation), vegetation factor (enhanced vegetation index), and bioclimatic factors (bio 
1-19). The bioclimatic factors were the following: Biol: Annual Mean Temperature; 
Bio2: Mean Diurnal Range; Bio3: Isothermality; Bio4: Temperature Seasonality; Bio5: 
Max Temperature of Warmest Month; Bio6: Minimum Temperature of Coldest Month; 
Bio7: Temperature Annual Range; Bio8: Mean Temperature of Wettest Quarter; Bio9: 
Mean Temperature of Driest Quarter; Biol0: Mean Temperature of Warmest Quarter; 
Bioll: Mean Temperature of Coldest Quarter; Biol2: Annual Precipitation; Bio13: 
Precipitation of the Wettest Month; Biol4: Precipitation of the Driest Month; Bio15: 
Precipitation Seasonality; Biol6: Precipitation of the Wettest Quarter; Biol7: 
Precipitation of the Driest Quarter; Bio18: Precipitation of the Warmest Quarter; and 
Biol9: Precipitation of the Coldest Quarter. Elevation data were obtained using the 
digital elevation model (DEM) of the Shuttle Radar Topography Mission (SRTM), 
which has a 90 m spatial resolution. The Enhanced Vegetation Index (EVI) map was 
obtained from Moderate-Resolution Imaging Spectroradiometer (MODIS) data 
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(https://modis.gsfc.nasa.gov/). Bioclimate factors were obtained from the WorldClim 
database (Fick and Hijmans, 2017; http://www.worldclim.org). All layers were 
resampled based on a bilinear interpolation in order to match the resolution. 


In order to prevent overfitting and identify the most appropriate variables, a correlation 
analysis was performed using Pearson's correlation coefficient (r) to examine the 
relationship between the 21 environmental variables. A Pearson's correlation value 
greater than 0.7 signifies a strong correlation between the two variables. All but one of 
any given environmental layers with a high correlation with other layers were removed. 


For future climate, bioclimatic variables using the MIROC-ESM climate model for the 
period 2041-2060 were downloaded from the WorldClim database 
(http://www.worldclim.org). Two representative concentration pathways (RCP), namely 
RCP 2.6 and RCP 8.5, were used as scenarios to predict future climate change. The 
RCP2.6 and RCP8.5 scenarios reflect low and high greenhouse gas emissions in the 
future, respectively. 


Species Distribution Modelling 


The species distribution models were created using MaxEnt v. 3.4.4 (Phillips, Anderson 
and Schapire, 2006). MaxEnt modeling was conducted by allocating 25% of the 
occurrence data for testing purposes and 75% of the occurrence data for training 
purposes. In order to identify the primary environmental factors that influence the 
distribution of P. stejnegeri habitat, a jackknife permutation analysis was performed to 
rank the environmental factors based on their percentage contribution. In order to assess 
the effectiveness of the MaxEnt model, the area under the receiver's operating curve 
(AUC) was employed as a measure of the model's predictive accuracy. A model's 
prediction accuracy increases as the Area Under the Curve (AUC) value increases. 
Evaluation of model performance were based on three specific criteria: poor, defined as 
an AUC value below 0.8; good, defined as an AUC value between 0.90 and 0.95; and 
excellent, defined as an AUC value above 0.95. Once the model computations were 
finished, the habitat distribution of P. stejnegeri was displayed on ArcMap. 


The continuous habitat suitability maps of P. stejnegeri for both current and future 
climate conditions (RCP 2.6 and RCP 8.5) were transformed into presence/absence 
binary maps using a 10 percentile training presence clog threshold. In order to evaluate 
the changes in the distribution of P. stejnegeri due to changing climate, binary maps 
representing the current climate conditions were overlaid onto the projected habitats and 
compared using the SDMtoolbox in ArcMap. Afterwards, the results were categorized 
into three distinct groups: (1) range expansion, (2) no occupancy (absent in both), (3) no 
change (present in both), and (4) range contraction. 


Results 

Variables Correlation and Selection 

The correlation analysis of all 21 predictor variables showed a high correlation for some 
variables, as shown in figure 2. The variables with the highest correlation are Bio6, Bio9, 


Biol0O, and Biol! in relation to Biol, which has a correlation coefficient above 90% 
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positive correlation. Nine variables were selected out of 21 variables after the variable 
selection process, namely Biol (Annual Mean Temperature), Bio2 (Mean Diurnal 
Range), Bio3 (Isothermality), Bio4 (Temperature Seasonality), Biol2 (Annual 
Precipitation), Biol5 (Precipitation Seasonality), Biol8 (Precipitation of Warmest 
Quarter), Biol9 (Precipitation of Coldest Quarter), and elevation. 
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Figure 2. Correlation analysis of 21 predictor variables. a.) Heat map showing the 
correlation between variables, b.) Dendrogram showing the grouping of highly 
correlated variables with correlation coefficient threshold of 0.7. 


Key Environmental Factors Influencing the Habitat Distribution 


The importance and contribution of each environmental variable to the prediction model 
were analyzed by the jackknife method. The results showed that the contributions of 
Biol (Annual Mean Temperature), Elevation, and Biol2 (Annual Precipitation) were 
24.39%, 18.65%, and 17.21%, respectively (Table 1). The environmental variable with 
the highest gain when used in isolation is Biol (Figure 3), which, therefore, appears to 
have the most useful information on the distribution of P. stejnegeri. The environmental 
variable that decreases the gain the most when it is omitted is Bio4, which appears to 
have the most information that is not present in the other variables. 


Table 1: Estimates of relative contributions of the environmental variables to the MaxEnt 
model 


Variable Percent contribution Permutation importance 

Biol 24.39 28.31 
Elevation 18.65 2.026 

Biol2 17.21 11.42 

Bio4 14.17 18.61 

Bio2 9.47 18.19 

Biol8 8.2286 19.5545 

Biol9 7.0256 0.1321 

Bio3 0.4861 0.1612 

Biol5 0.3655 1.5926 
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Figure 3: Jackknife test of variable importance on P. stejnegeri distribution 
Modelled Distribution of P. stejnegeri under Current Conditions 


The Maxent model predicts the distribution of P. stejnegeri across the entire island, as 
shown in figure 4. The species’ current modeled distribution closely matches the region's 
published distribution data, including the range map produced by the IUCN 
(International Union for Conservation of Nature, 2024). The presence/absence binary 
map showed that the entire area where P. stejnegeri is predicted to be present is 
38,550.18 km?. The model AUC value of 0.917 indicates that the model's performance 
is superior to random prediction, hence confirming its validity. 
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Figure 4: Species distribution model of P. stejnegeri under current conditions. 
a.) MaxEnt model of P. stejnegeri showing probability of occurrence; 
b.) Presence/absence map of P. stejnegeri; c.) Receiver operating characteristic (ROC) 
curve for evaluating the model performance of MaxEnt model. 


Neil Jun Sala Lobite 128 


Grassroots Journal of Natural Resources, Vol.7, No.1 (April 2024), p.123-137 | ISSN 2581-6853 | CODEN GJNRA9 
Doi: https://doi.org/10.33002/nr258 1 .6853.070107 Open Access 


Neil Jun Sala Lobite 


Modelled Distribution of P. stejnegeri under Future Climate Scenarios 


Using the MIROC-ES model, the distribution of P. stejnegeri is projected under future 
climatic scenarios. Both models under RCP 2.6 and RCP 8.5 indicated a decrease in the 
species’ range, as shown in figure 5, using the 10th percentile threshold cutoff. The 
predicted extent of P. stejnegeri's present under RCP 2.6 is 17,890.26 km’, but in RCP 
8.5 itis 8,452.72 km”. The findings indicated a decrease in the extent of the area occupied 
under projected future climate conditions. 
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Figure 5: Species distribution model of P. stejnegeri under future conditions (2050). a.) 
MaxEnt model showing probability of occurrence under RCP 2.6; b.) MaxEnt model 
showing probability of occurrence under RCP 8.5; c.) Presence/absence map under 
RCP 2.6; d.) Presence/absence map under RCP 8.5. 


Figure 6 illustrates the geographical areas where the distribution of P. stejnegeri will 
expand and contract in response to changing climate conditions, as determined by 
different RCPs. According to the RCP 2.6 climatic scenario for the year 2050, the 
suitability for P. stejnegeri would decrease, resulting in a total range contraction area of 
23,164.77 km? (Table 2). The majority of the contractions are located in the northern 
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region of Mindanao. On the other hand, under the extreme climate scenario (RCP 8.5), 
the area of suitable habitat is predicted to be greatly reduced, with a total contraction 
range of 30,200.16 km’, or 31.72% contracted in reference to the current model. The 
areas where no change occurs, 1.e., species that are predicted to be present in both current 
and future climatic scenarios, are primarily concentrated in the high-elevation parts of 
Mindanao. 
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Figure 6: Distribution changes of P. stejnegeri from the present climate to the future 
year (2050) under climate scenarios: (a) RCP 2.6; (b) RCP 8.5. 


Table 2: Relative changes in P. stejnegeri distribution area under climate scenarios RCP 
2.6 and RCP 8.5 for the year 2050. 


Distribution Total area of change (km’) Percentage of change 
Change RCP 2.6 RCP 8.5 RCP 2.6 RCP 8.5 
Range 2,519.68 101.99 2.65 0.11 
expansion 
No occupancy 54,116.55 56,534.24 56.84 59.38 
No change 15,400.04 83,64.65 16.18 8.79 
Range 23,164.77 30,200.16 24.33 31.72 
contraction 
Discussion 


The present distribution modeling of P. stejnegeri is crucial for understanding the 
environmental factors that influence its distribution. The MaxEnt model precisely 
delineates the habitat range of P. stejnegeri, exhibiting a close correspondence with 
established distributions and conservation evaluations, such as those undertaken by the 
IUCN. The high AUC score validates the model's accuracy in identifying the significant 
environmental factors that influence the distribution of P. stejnegeri. The importance of 
the average annual temperature, altitude, and annual rainfall in the model emphasizes 
the species' reliance on these biological parameters. This is consistent with overarching 
ecological theories that credit temperature and precipitation as the primary factors that 
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influence patterns of species distribution (Bickford et al., 2010; Guisan and Thuiller, 
2005; Mi eft al., 2022). The mean annual temperature significantly influences the 
distribution of amphibians, as these organisms, being ectothermic, depend on external 
heat sources to regulate their body temperature. The metabolic rates, reproductive cycles, 
and survival of frogs are intricately linked to the prevailing ambient temperature 
conditions (Rollins-Smith and Le Sage, 2023). Wake and Vredenburg (2008) have 
provided evidence of the significant impact that changes in temperature patterns can 
have on amphibian populations. The vulnerability of P. stejnegeri to the annual mean 
temperature underscores the species' sensitivity to climate change, hence reflecting 
broader concerns regarding the impact of global warming on amphibian populations. 
Additionally, annual precipitation is another important factor that indicates the essential 
requirement for moisture by amphibians. Due to their permeable skin, these organisms 
rely on a moist environment to facilitate gas exchange (Harvey Pough, 2007). This 
unique physiological trait makes them highly susceptible to fluctuations in moisture 
levels within their habitat (Harvey Pough, 2007; Lips et al., 2003; Stuart et al., 2004). 
The study conducted by Mi ef al. (2022) clearly illustrates the correlation between 
precipitation patterns and the availability of suitable habitats for amphibians, showing 
that the average annual precipitation had the greatest impact on the spatial distribution 
pattern. Based on the model, the reliance of P. stejnegeri on the annual precipitation 
factor further demonstrates the complex interplay between amphibians and their 
hydrological habitats. Lastly, elevation was found to be a significant factor for P. 
stejnegeri's distribution. It has been proposed that elevation affects the local weather 
conditions, which, in turn, affects amphibian population diversity and dispersal. McCain 
and Grytnes (2010) demonstrate how elevation gradients influence biodiversity patterns, 
such as those of amphibians, by generating diverse ecological niches at different 
altitudes. The importance of elevation in P. stejnegeri's distribution model emphasizes 
how changes in altitude can define the boundaries of acceptable habitats for amphibians, 
affecting their distributional range and population dynamics (Fu et al., 2006; Supsup et 
al., 2022). Overall, the model's ability to accurately predict the current distribution of P. 
stejnegeri suggests that the species is highly influenced by a range of climatic and 
topographic parameters, specifically temperature, precipitation, and elevation. This 
highlights the crucial significance of these environmental conditions in the conservation 
of this species. 


According to future climate projections, it is predicted that the range of P. stejnegeri 
would see substantial changes, with a large reduction in appropriate habitats, especially 
under the high-emission scenario (RCP 8.5). This estimate is consistent with the global 
consensus on the anticipated effects of climate change on biodiversity suggesting that 
several species would have reduced habitats because changing climatic conditions 
(Bellard et al., 2012; Blaustein et al., 2010; Chen et al., 2011). The model forecasts a 
significant decrease in the spatial range of P. stejnegeri under both RCP 2.6 and RCP 
8.5 scenarios, highlighting the severe threats that climate change poses to its habitat 
suitability. The range contractions in the northern part of Mindanao are particularly 
noteworthy, suggesting that these areas may no longer be suitable for the species due to 
shifts in climatic conditions. The anticipated decrease in habitat in Mindanao could have 
significant implications for the population of P. stejnegeri, given that this species is 
endemic to the island. The projected changes in the distribution of P. stejnegeri 
emphasize the broader challenges that endemic species in biodiversity hotspots are 
facing because of climate change (Trew and Maclean, 2021). In light of these concerns, 
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potential strategies could involve habitat restoration, the creation of biological corridors, 
and the safeguarding of climate refuges that may remain suitable amidst changing 
conditions (Rahman et al., 2019; IUCN SSC Amphibian Specialist Group, 2020). The 
establishment of corridors will be useful in linking habitats and refuges of P. stejnegeri 
that are expected to remain viable under future climate scenarios. This will facilitate the 
movement of the species to more suitable environments and is crucial for its adaptation 
to environmental changes (Beier, 2012). 


Conclusion 


This study models the habitat suitability and distribution of P. stejnegeri under both 
current and projected future climate scenarios using MaxEnt. The study determined that 
annual mean temperature, elevation, and yearly precipitation are the main factors 
affecting P. stejnegeri's distribution. Future projections demonstrate the vulnerability of 
P. stejnegeri to climate change, highlighting a substantial decrease in its geographical 
range under both mild (RCP 2.6) and severe (RCP 8.5) future climate conditions. The 
reduction in range, specifically in the high-emission RCP 8.5 scenario, emphasizes the 
pressing risk that climate change poses to the species. Given the species’ sensitivity to 
environmental factors such as temperature and precipitation, together with the projected 
decrease in geographical range due to climate change, it is important for conservation 
efforts to prioritize the adaptability and quality of P. stejnegeri's habitats. 
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